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In this article we discuss the interference patterns which appear in the local particle- 
density and particle-current distributions of two infinite 2-D lattice systems of free fermions 
which are allowed two communicate via two tunneling junction points. The two fermion 
reservoirs are initially in different invariant states, and the particle-density and particle- 
current distributions are calculated in the final stationary state. The dependence of the 
interference pattern on the interaction details is discussed in several examples. 

1 Introduction 

The paper is concerned with a simple model of electron tunneling between infinite reservoirs in 
0\ ', different equilibrium states. We work in a tight binding approximation, i.e. the electrons live 
on lattices and the one-particle kinetic energy is given by appropriate lattice Laplace operators. 
Thereby, the interactions are neglected, so that the Hamilton operators of the isolated reservoirs 
are bilinear in the electron creation/annihilation operators. At time t = a coupling between 
O . the reservoirs is switched on, likewise bilinear, which allows direct hopping of the particles 
between certain finite subsets of sites in different reservoirs. Asymptotically in time, the system 
approaches a stationary state, in which permanent currents of particles and energy are present. 
^ ' For such models, the approach to stationarity and the structure of the stationary state is well 
■ understood at the mathematical level of rigor [1], |2]- The stationary state is a quasi- free state 
[3], which is completely characterized in terms of the spectral properties of the one-particle 
Hamilton operators of the free and coupled reservoirs. Explicit formulas for the permanent 
particle and energy currents between systems are known (the Landauer-Biittiker formulas). 
However, the state provides more than such global information about the system: it allows 
obtaining the expectation values of all the local observables. 

Our aim here is to perform such a detailed study of the stationary state, by calculating the 
expectations of various local observables (namely, the local particle density and the particle 
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currents) for a particular geometric arrangement and tunneling Hamiltonian. More precisely, 
we take the reservoirs to be free fermion systems living on regular two-dimensional lattices, 
with direct tunneling between (one- and/or two-site) finite subsets of the latter. Such choice of 
the reservoirs is by no means unrealistic, e.g. electrons in the surface states of the close-packed 
metal surfaces are two-dimensional, nearly free, electron gases. Moreover, the local density 
of electrons is experimentally accessible by scanning tunneling microscopy, by a technique 
elaborated in the nineties by Crommie et al |4j. Standing density waves were observed on the 
Cu(lll) surface and their origin explained by the scattering on point impurities located on the 
surface |1],[5]. Such measurements have since been performed on different materials and defect 
configurations (see [6] for a recent experiment), and various applications proposed (e.g. [7]). 
Though our calculations refer to a different physical situation, we hope that they are not only 
of theoretical interest and that the interference effects exhibited by us can be experimentally 
observed. 

The paper is organized as follows. The precise description of the model and its general 
properties are done in Section [2l In order to establish notation and for reader's convenience, 
the approach to, and structure of, the stationary state are presented in Subsection 12. ![ following 
[Ij, [2j. The main ingredient is the perturbation analysis of the coupled Hamiltonian. As we 
want to demonstrate the appearance of quantum interference effects in the stationary state, 
we derive the formulas for the particle density profile (Subsection 12. 2p and the local particle 
currents (Subsection l2.3l) in the stationary state. Section |3]is devoted to numerical illustrations. 
Here, we specify the reservoirs to be two-dimensional free fermion lattice gases, with direct 
tunneling between (one- and/or two-site) finite subsets of the two square lattices. Subsection 
13.11 contains the numerical study of the limit values on the real axis of the Green's function 
of the two-dimensional lattice Laplacian, which enter the expression of the stationary state. 
Subsections 13. 2[ 13.31 are devoted to numerical calculations of the charge- and of the current- 
density profiles, respectively, for the particular geometry described above. Classically, the total 
current should be a sum of the currents established through each individual channel when the 
other are closed. The computations show that this is no longer true in the quantum case. 
Moreover, the calculated charge-density profiles in the reservoirs exhibit interference fringes, 
depending on the tunneling constants for every pair of coupled sites. Section H] is devoted to a 
discussion of the results. 

2 The stationary state 

2.1 The model and the approach to stationarity 

We consider two lattice free electron gases on Lj = Z'^', i = 1,2 (where Z denotes the integers), 
coupled via a tunneling junction of finite support. A one-electron wave function is a vector 
f e n = fiU U L2) = /2(Li) © l\h2), i.e. / = {/^; x G Li U L2}, such that := 
X] < Let a*, a J. be the operator of creation, respectively annihilation, of an 

3:GLiUL2 

electron at the site x acting in the anti-symmetric Fock space J-" over H. Also, for / G 
we denote a*(/) = J2x f^'^x ^(/) — 'I2xf^^x- ^ one-particle operator c, given by the 
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matrix {cx,y) {x,y e Li UL2), its second-quantization, acting in J^, is defined as 

a;,2/eLiUL2 

Putting the electron mass equal to 1, the one-particle Hamiltonian corresponding to 
uncoupled reservoirs, which represents the kinetic energy operator, is given by 

2 2 

{h'f). = Y.{h'j%^-\Y.XK{^) E (A-/-)' (2-1) 

•i=l i=l ye'Li,\y—x\=l 

where Xa{^) denotes the characteristic function of the set A. The evolution in time (in the 
Heisenberg picture) of a^{f) writes as 

a?(a«(/)) := e'*'^r('^°)a«(/)e-'*^^('^°) = a«(e^'^°*/)- 

The generalized eigenfunctions (t){k) of are plane waves in each of the two reservoirs, 
corresponding to given momentum A; e Ti U T2, where Tj = (— tt, tt]*^*: 

ForA;eTi, (t){k):, = xk{x){'^t')~''''^ ^^V {ikx), (z = l,2). (2.2) 

One has hP(t){k) — uj{k)(f){k), where the dispersion law is 

2 di 

^{k)^J2 X-^i (^) E 2 sin' kj2. (2.3) 

1=1 a=l 

At time t < 0, the two (uncoupled) reservoirs are supposed to be in different equilibria, 
defined by inverse temperatures /3j and chemical potentials /ij, i = 1,2. Let < • >o denote this 
constrained equilibrium state. The density operator of the state < • >o is 

2 2 
{p\,y:=<a:ay>o=Y^ J dkmMWyfPi,,Mk)) = Yl / dkE{k),,yf^„,Mk)): (2.4) 

where fj3 nix) — {1 + e^(^~''))~^ is the Fermi function, and where we defined the matrix E{k) 

by ' 

E{k)^,y = ^{k),^^{k)y, x,yGLiUL2. (2.5) 
The higher order correlations in the state < • >o are calculated in terms of these as: 

< a*{fi)...a*{fm}a{gi)...a{gn) >o= Sm,ndetM, M^j =< a*{fi)a{gj) >o= {gj,p^fi)H- (2-6) 

A state with the latter property is called quasi-free and gauge-invariant. 

Starting from t = 0, a tunneling junction is opened, connecting the finite subsets Si G Lj, 
i = 1,2. The future time evolution is controlled by the perturbed one-particle Hamiltonian 
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h = + v, where v = t + t* with r : P{h2) — t- /^(Li) describing the tunnehng of particles from 
the second to the first reservoir, i.e. 

irf). = XsAx)J2t,,yf{y), (2.7) 
y&S2 

and T* : /^(Li) — )■ /^(L2) the tunnehng in the opposite direction {t*)x,y = ty^x- The Heisenberg 
evolution of a^{f) for t > writes 

a,(a»(/)) := ei*dr(/^«+.)^fl(j)e-iMr(/.o+.) _ aS(ei('^"+-)*/). 

The initial state < ■ >o is no longer invariant under the new evolution. At time t > 0, the 
state will be the gauge- invariant, quasi-free state < ■ >t of density operator 



-\th° AthT.* J) -\th° Ath . 



(2.8) 



here, the factors e'*'^" are added for convenience (they commute with hence do not con- 
tribute), as they show that, in far future, the state is given by the Moller operators in the 
scattering theory for the perturbed evolution. More precisely, we shall use the following result, 
which is an easy consequence of the general results of [T], [2]: 

Proposition 2.1 In the framework defined above, 

1. the spectrum of h consists of absolutely continuous spectrum, equal to the spectrum of , 
o"ac(^) = 0"(/i°) = Uj=i.2[0, 2(ij] and of a finite set of eigenvalues, o"p(/i), outside (y{h^); 
thereby, the wave operators 

W±:=s- lim e'*^-'*'''' (2.9) 

exist and establish the unitary equivalence of with the absolutely continuous part hac of 
h: hacW± = W^h^; 

2. the states < ■ >t approach in mean, as t ^ oo, a stationary state, i.e. the following limits 
exist: ^ 

lim (1/T) [ < A >i dt =:< A >stat, (2.10) 
for all polynomials A in a^{f); 

3. < ■ >stat is the gauge-invariant, quasi-free state of density operator 

p+ = W.p'>Wl+ J2 ^aP°Pa, (2.11) 

Ag(Tp(h) 

where Px is the projection onto the eigenspace of h corresponding to the eigenvalue X. 
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The rest of this section is devoted to deriving computable expressions for the eigenvectors 
and of the Moller operators. 

We start with writing the resolvent operators r°(2;) = {h^ — z)^^ and r{z) = {h — z)^^. 
Clearly, r^{z) is expressed in terms of the generalized eigenfunctions (12. 2p of as: 

AzUy = [ 4^dfc = / ^^de, (2.12) 
J uj[k)-z J e-z 

T1UT2 

where E{k) is defined in Eq. (l2.5p and P{e) is the integral of E{k) over the energy shell Si(e) U 
E2(e), with 

Si(e) = {ke Ti- uik) = e}. 

Both E(k) and -P(e) have block diagonal structure with respect to the reservoirs. In particular, 
P(e) = -Pi(e) © -P2(e) where, denoting d(jj(fc) the element of surface area on Sj(e), we have, 
explicitly, 

ne).,, = X^XL.(x-y)P.(e).,, = X;XL,(a:-y) J (2-13) 

S,(e) 

2 

Thereby, r°{z)^^y = Yl Xui^ - y)9i{z^ V - a^), where gi{z, ■) : Lj -> C is defined as 

i=l 

9.{z.x)= [ ^^de. (2.14) 

J € Z 


To calculate r{z), one has to solve for / the equation (/i° — z + v)f = g. li z E C\ cr{h^), this is 

equivalent to f + r^{z)vf = r^{z)g. Applying v to the latter equation, one gets {v + vr^ {z)v) f = 
vr'^{z)g. Let 11 be the projection on the range of v: IVH = (kerf)-*-. There are two possibilities: 

• The restriction of v + vr^{z)v to HH is invertible. Then, 11/ = {v + vr^{z)v)^^vr^{z)g. 
Once n/ is known, the whole vector / is / = —r^{z)v{Ilf) + r^{z)g. For such values of 
z, h — z has a bounded inverse r{z): 

r{z) = r°(;z) - r°{z)Q{z)r^{z), where Q{z) = v{v + vr'^{z)v)~^v; (2.15) 



• There exists a non-zero vector ip G TYH, such that {v + vr^{z)v)il) = 0. Such values of z 
are eigenvalues of h with eigenvector / = —r^{z)vil!. Conversely, if / is an eigenvector 
of h for the eigenvalue z, then 11/ 7^ and / + r^{z)vf = 0. Applying v, it follows that 
{v + vr^{z)v)f = {v + vr^{z)v){Ilf) = 0. Therefore, (Jp{h) equals the set of z with this 
property and there is a one-to-one correspondence between the eigenspace corresponding 
to z and UH fl ker [v + vr'^{z)v). 
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We consider next the Moller operator entering Eq. (12. lip . As a prerequisite, we remark 
that the functions gi{-; x) defined in Eq. fl2.14p are analytic in the complex plane cut along the 
segment [0, 2di\ and have finite boundary values, i.e. the limits 

\im gi{e ± ie, x) =: gi^±{e;x) (2.16) 

exist for all e G (0, 2di) and x G L^. As a consequence, and because v has finite range IIH, the 
limit values 

Q±(e) := limQ(e ± ie) 

exist. Clearly, Q±{e)x,y = unless x,y & 5*1 US'2, so we shall view Q±{e) as a matrix indexed by 
5*1 U 5*2, and denote Q±{e)x'y^ its submatrices corresponding to x E Si and y G 5*^ = 1, 2). 
Applying W^, Eq. (12.91) . to a generalized eigenfunction of h^, one has: 

[W.(j){k)]x = limj^ e {6x, e-"(''-^('=)-^^V('^)) dt 
= lim {6x, {—ier{u{k) + ie)(j){k)) 

= - [ri{u{k))Q+{uj{k))m]x ^^-^^^ 

= {6x,[I-rl{uj{k))Q4oo{k)Mk)), 

where we used Eq. (l2.15l) for r{z) and the obvious relation {—ier^{u{k') + ie)(j){k') = (f){k'). 
Likewise, 

[wimi = (^x., [/ - Q-i^ik)y_iujik))]m)- (2.18) 



We shall demonstrate the interference effects on two properties of the stationary state Ostat? 
the permanent currents and the particle density. In the calculations below, we make the as- 
sumption that tx^y (hence, all hx^y) are real numbers (physically, complex h^^y are related to the 
presence of magnetic fields). 



2.2 The permanent currents in the stationary state 

Let Ni^x} = dx^x = dr{xx) be the particle number operator at the site x G Li U L2, where Xx 
is the projection onto the vector 6x, i.e. its matrix is {Xx)y,z = ^x,y^x,z- If the evolution is given 
by the one-particle Hamiltonian h, the fiux of particles from the site x is defined as 

d 

I{x} ■= --^aMT{xx)\t=o = -i[dr(/i),dr(xx)] = dr([-i/i, xj)- 

In the stationary state, as < at{N{x}) >stat is independent of time, we have the "continuity 
equation", saying that the total current fiowing from the site x vanishes: 

J{x} ■=< dT{[-ih,Xx]) >stat= 0. 
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The l.h.s. of this equation equals ^ 253< hy^;j.a*ax >stat- We define the current along the 

?/eLiUL2 

oriented bond {x, y} in the stationary state as: 

jx,y = 255< hy^^aya^ >stat = '^hy^^'^^d^, p+5y). (2.19) 

Clearly, j^^y = —jy,x and jx,y = unless h^^y 7^ 0, what happens, under our assumptions, only if 
either h^y ^ 0, implying that x,y are nearest neighbors in one of the two lattices, or Vx,y 7^ 0, 
implying that x E Si and y & S2, or x E S2 and y E Si. 

Note that, in view of the fact that all hy^^ are real, the eigenvectors of h can be chosen to 
have real components. It follows that the point spectrum (i.e. the second term in Eq. f|2.1ip ) 
does not contribute to j^^y, and Eq. f|2.19p writes 

jx,y = 2hx,y'^{W16x,p''W16y)=2hx,y ^ p^k)^ {{W 16 ^ , E {k)W 16y)}dk . (2.20) 

T1UT2 

Let A C Li UL2, and Xa = Xx be the multiplication by the characteristic function of A. 

xeA 

Then, dr{xA) = X] (^x^^x ='■ N\ is the number of particles in A. The fiux of particles out of A 
is 

Ja := ^ -\hx-y{alay - ala^), 

(x,y)edA 

where the /i-boundary of A, dA, is the set of ordered pairs (x, y) of sites x, ?/ G Li U L2 such 
that hx,y 7^ 0, x e A, y G Li U L2 \ A. 

If A is finite, then < at{N/<^) >stat is finite and independent of t, therefore Ja :=< I\ >stat= 

E jx,y = 0. 
(x,y)edA 

If A is infinite, then < A^a >o is infinite, whenever the particle-density is non-zero in the 
initial state, i.e. /ij > in Eq. ( 12.41) . The same holds for the evolved states < • >t and their 
hmit < ■ >stat- 

We consider below the case A = Li, whereby we suppose pi > 0. Then, A^Li repre- 
sents the number of particles in the first reservoir, which is a conserved quantity for i.e. 
[A^Li,dr(/i°)] = 0. However, in this case, dhi = {{x,y}; x G 5*1, y G 5'2} is finite, so the 
stationary current from the first reservoir Jl^ =< Jl^ >stat makes sense and equals 

x£Si,y£S2 

whereby h^^y = v^^y. It is worth mentioning that, for any finite sets Ai,A2, such that 5*1 C 
Ai C Li and 6*2 C A2 C L2, the following equalities hold: 

-'li\Ai = '^LiUA2 = 

showing that the same current traverses any finite contour surrounding 5*1 in Li, or surrounding 
5*2 in L2. 
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We shall calculate below j^ y, Eq. (l2.2Up . for x G Si, y G S2, using Eq. (I2.18p . One has 

{W*_S^,Eik)WlSy) = iS,,[I-rlQ+]Eik)[I-Q_r^]Sy), 

where Q±,r^ are calculated at e = uj{k). Using that both r^(A;) and E{k) have block-diagonal 
structure with respect to the reservoirs, one obtains 







)(1.2)^0 



+2v,,y J dkpl'^ -rl^Q^i''^ E,{k) + rl^Q'l^-'^ E2{k)Q'^^-'Wl 



n(i.2) 



)(2,2) n 



(2.22) 



Taking advantage of the fact that all factors, but Ei{k), under the integral signs depend solely 
on e = u}{k), one can perform the integrals over the energy shells Sj(e) as in Eq. (12. 13^ : 



/x,y 



2di 



(1,1)/ 



x,y 



(2.23) 



2d2 



+2vy, J depO(e)^> -rO + (e)QV'''^(e)P2(e) + rO + (e)gV'"^(e)P2(e)QL''^^(e)rO _(e 



(1,2), 



^(2,2)/ ^ 



x,y 



The summation over x,y in Eq. (l2.23p allows a significant simplification of the expression 
for the total current: 



Jl, = 2QJ dep;(e)tri ( -Pi(e)gL'''^(e)rO_(e)t; + r?.+ (e)g;'')(e)Pi(e)g^'^(e)rO „(e)t; 
^ 

+25>'/ dep0(e)tr2 (-t;rO + (e)Qi_''')(e)P2(e) + t;rO + (e)gi!''^(e)P2(e)QL'''^(e)r20_(e) 
^ 

2(ii . . 

= 253/ dep?(e)tri (Pi(e)QL'''^(e) - Piie)Q^l^'\e)rl4e)Q^^''\e)) 
^ ^ 

+25>/'depO(e)tr2 (Qf'\e)P2ie) - Q^''''\e)rl4e)Qf^'\e)P2ie 



where we applied the identity Q{z)r^{z)v = vr'^{z)Q{z) = v — Q{z) and the cyclic invariance 



of the trace. Now, 2'^tTA = -itr(A - A*) and Q^l''^ - Q"^^'' = E,- <5+^^(^i,+ - 



)(*,i)/ 



i(i,«) 



E,Qf'^2i7rP,QL^'^so 



(i,i) 



2 mill (di ,(i2) 



Jli = 27r J de(p?(e) - p°(e))tri (Pi(e)gL''^)(e)P2(e)Qj'^)(e)) . (2.24) 
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2.3 The charge- density profile in the stationary state 

The charge density at the site x G L1UL2 in the stationary state is given, according to Eg. (12. lip , 
by: 

d{x) :=< a*.a^ >stat= (5x,P+5x) = dp{x) + d^c{x), (2.25) 
where we separated the contributions of the point and absolutely continuous spectra: 

AGo-p(/i) j=l,2fceTi Aecrp(h) 

(2.26) 

i=l,2A;eTi 

Clearly, E dpi^) < 00, meaning that the number of particles accommodated on the eigen- 

states of h is finite and located near SiU 32- 

We calculate dg^d^) for a; € L2. Remembering Eq. (l2.17p and the definition (12.50 . we have 

I [W.mi r = {^x, [I - rl{co{k))Q+iu{k))]Eik)[I - Q^{co{k)y_{co{k))]6,) . (2.27) 

Again, the integration over the energy shell uj{k) = e can be performed with Eq. (I2.13p . so, 
finally, 

d.c{x) = E / de/ft,M.(e) (5., [/ - r° (e)g+(e)]P,(e)[/ - Q -{e)r'_{e)%) . (2.28) 



i=l,2 







Here, one should take into account the block structure of r° and the fact that Pi{e)y^z = 
unless both y,z & Lj. For instance, if x G L2, the term i = 1 in Eq. (l2.28p simplifies to 

S'\x) = J def^,,,,{e) J2 92Ae;s2-x)[Q^^''\e)P^{e)Q^}'\e)]s,,s',92,-{e;x-s'2), (2.29) 

while the term i = 2 has a more complicated structure, due to the superposition of the incident 
and the refiected waves: 

2d2 

S'\x)= / de/;32,M2(e){P2(e).,. 




^72,+ (e;x-S2) (Qf'\e)P2{e)) + (P2{e)Q^'''\e)) ^72,-(e;s2- 



X 



- E 

s2eS2 

+ E 92Ae; X - S2)[g?'^(e)P2(e)gL''')(e)],,,,^72,-(e; 4 - x)}. 

(2.30) 

The following remark is in order. If the dimension ^2 > 2, the Green function 5'2,±(e; x) — )■ 

2d2 

as X — i- 00. Hence, as x —i- 00 in L2, d{x) approaches the constant P2 = J deff3,^^^2i^)P2{G)x,x, 



which is the initial equilibrium density of the second reservoir. Indeed, turning on a tunneling 
junction should not affect the intensive parameters defining the reservoirs, this is one of the 
reasons why we restrict to dimensions of the lattices larger than 1. 
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3 Numerical computations 



In the computations we particularize to two-dimensional homogeneous reservoirs, i.e. we take 

Li, L2 both equal to T? , and h\, h\ both equal to -(1/2) A acting in t^il?) , see Eq.([2ll]). 

2 

Hence, Ti = T2 = T = (— tt, tt)^ and the dispersion law of is u{ki, = X] 2 sin^ {ka/2). 

a=l 

3.1 Evaluation of the Green's function of the 2D Laplacian 

We consider here the function g{z; x) = [(—(1/2) A — z)^^]^^. It is defined, for z in the complex 
plane cut along the real segment [0, 4] and for x = (xi, X2) € Z^, by the formula (I2.12p : 

/Q-i(kiXi+k2X2) 
dkidk2^j — . (3.1) 

T 

Remark 3.1 As done in Eq A2.1^) . we can integrate first over the energy shells. For e < 2, the 
energy shell S(e) = {k = {ki, ^2) ^ T; u{k) = e} is a closed curve, which can be parametrized by 
one angle 9 G [0,27r). Indeed, for all e G (0,2) and any 9, there exists one point k{e,9) G S(e) 
along that direction. Denoting K{e, 9) the norm of this k, we have 

ki{e,9) = K{e,9) cos9 
k2{e,9) = K{e,9)sm9 ' 

with K{e, 9) the unique solution of the equation u{K cos 9, K sin 9) = e. The change of variables 
{ki,k2) — )■ {e,9) is regular from \ki\ + 1^2! < to (0,2) x [0, 27r) and its Jacobian equals 

J(e, 9) = K{e, 9)deK{e, 9). (3.2) 

Hence, if e < 2, and denoting (j) the polar angle of x = {xi,X2), we have 

P{e)o,x = j d9J{e, 9) exp [-i\x\K{e, 9) cos (0 - 9)]. (3.3) 



For e > 2, S(e) consists of four disconnected parts, which can be put together by changes of 
variables of the form /c^ = fc^ ± vr into the closed curve corresponding to S(4 — e) (see Figure 
H]), which can be parametrized by an angle 9 as before. 

One can easily see that 

g{z; xi,X2) = g{z; X2, xi) = g{z; |a;2|) = (-l)^+'''+^"5'(4 - z;xi, X2) = g{z; Xi, X2). (3.4) 

Hence, the values of g{z] x) are real for real z outside the segment [0, 4]. We plotted in Figure 
[21 as an example, '^g{z; 0, 0), which is strictly positive on the upper part of the cut. Also, g is 
determined by the values g{z; m, n) with < m < n and ^z < 2). 
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Figure 1: Energy levels S(e) 



Using the identity = i e-'^*dt, if < 0, g can be expressed, for z in the upper 
complex half-plane, as an integral of a product of two Bessel functions of integer index, 

^(^;m,n)=i™+"+i / dte-'*^^-) J„(t) J„(t), (3.5) 
Jo 

where we used the representation 



(-i)" r 

Jn(t) = / exp(— in/c + it COS A;)dA;. 

2vr 



The integral (13. 5p is related to a regularized hypergeometric function [8j : 



^ \ m+n+l ^ 



g{z;m,n) = - 2)"^+"+i ("^ + ^)' " 4^3'^(ai, ^2, "4; &2, ^s; -^j^^), (3.6) 

where 

1+m+n 2+m+n , ^ , ^ , ^ 

oi = a2 = ^ , as = a4 = ^ , Oi = 1 + m, 62 = 1 + ^, 03 = 1 + + 

4F3^^(...; ...;«) is analytic in the complex-M-plane cut along [1, 00). Note that the limit values 
g{e + iO; ...) are expressed in terms of 4^F^'^^{...; j^Zq^ — iO) if < e < 2. 

For m, n not too large (say, m + n < 5), Eq. ( 13.61) (including the limit values g{e + iO; m, n) 
for < e < 4) can be safely computed using the built-in functions of Mathematica. For larger 
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Figure 2: Graph of '^g{z; 0, 0) in the complex 2;-plane cut along [0, 4] 

values oi m + n, an asymptotic study of g has been performed, based on the steepest descents 
method [9], [10], the output of which is described below: 



For z outside [0, 4], g{z;x) is the Fourier transform of a C°°, periodic function of ki,k2, 
hence it decays faster than any inverse power of the distance 

When z = e + ie, < e < 4, the denominator in Eq. (13.11) develops a singularity in the 
integration domain when e \ 0. As remarked below Eq. (13. 4p . it is sufficient to consider e < 2 
and X = {m,n), m,n> 0. Let x{t) be a C°°-function, which equals 1 if |t — e| < a and vanishes 
if |t — e| > b, where 0<a<6<2 — e. Then, the asymptotic series of g{e + iO; x) coincides 
with that of 

,-2 Z"^,^ x(2-e-cosA:i-cosfc2)e-'(^^"+^^") 
,,(e + i0;x) = (2vr) J dk.dk, 2 - e - iO - cos fc, - cos fc^ " ^'"'^ 

T 

Therefore, remembering Eq. (13.31) . 

g^{e + iO; x) = j dt^^^^P{t)o,.. (3.8) 

The phase $(e, 0;^) = —K{e,6) cos {(p — 6) in Eq. (13. 3p (where (j) = arctan(n/m) is 
the direction of a; = (m,n)) has two stationary points with respect to 6, i.e. the equation 
dg^{e, 0; ^) = has two solutions 6s{e, 0) G [0, 7r/2] and 6s{e, (p) + vr, which control the asymp- 
totic behavior. Thereby, 

d^s^e,^;9)\e, = J{e,e,fij{e,(P) > 0, (3.9) 
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where 



ip{e, 0) = cos {Kg cos 6s) sin {Kg sin Og) sin + cos {Kg sin 6's) sin {Kg cos 6's) cos 0, (3.10) 

Here we denoted, for shortness, Kg = K{e, 6g{e, 0)), 6g = 6g{e, 0). Likewise, (9|$(e, 0; 0)|£i^_|_7r < 
0. Also, de^{e, 0; ^,) = -deK{e, Og) cos (0 - Og) < 0, as deK > and - 0, e [-7r/2, 7r/2]. 

A straightforward apphcation of Lemma 7.2.5 of [10] shows that the leading term of the 
|x| — )■ cxD asymptotical series of Eq. (13. 8p equals 

^(e + iO;x) / _ ^i(|x|j^,cos(0-e,)+V4) ^g^^^^ 



\x 



V27rV'(e,0) 



3.2 The density profile 







We compute here the density distribution (i.e. the number of band electrons per lattice site) 
in the second reservoir L2 in the stationary state. Thereby, we take the initial equilibria to 
correspond to zero temperature in both reservoirs, and to chemical potentials /ii in Li, and 
/i2 < /xi in L2. Hence, the initial (equilibrium) densities in the reservoirs take the values: 

pS = / ^(e)o,ode. (3.12) 

We choose for exemplification fii = 1.4, fi2 = 0.3, what corresponds to piq = 0.2804, peq^ = 
0.0492. 

Moreover, Si = {s*'^, s*'^}, {i = 1, 2) consist of two points at distance di apart along an axis 
of Lj, e.g. s*'^ = (0,0) and s*'^ = {di,0). The "contacts" between reservoirs is done by direct 
tunneling between the two pairs of points, {■s^'-', s^'-'} with tunneling constants tj, (j = 1,2). 
Thereby, we fix d2 = 20 and ^2 = 1, what proves to be good for a nice visualization of the 
interference phenomena, while di and ti are left as parameters. 

Under the assumptions made above, the transmitted and reflected densities in the second 
reservoir, Eqs. (12.291) and (12. 30 p . acquire the form (for simplicity, we do not exhibit the depen- 
dence on di and ti): 

ij'i 

rf«(x) = j 5^'\e;x)de {i = 1,2), x E L2, (3.13) 


where the functions under the integral sign have the following structure: 
i = 1 (transmitted) 

6^'\e;x) = (\/_(e,x),mtr(e)y_(e,a;))c2, (3.14) 

where 

mtr(e) = gj''^(e) ■ {P{e)\s,.s,) ■ Q^-''\e) > 0, (3.15) 
with P{e)\sixSi denoting the 2x2 matrix {P{e)g-^^g' , Si,s'i G Si} and V_(e,x) denoting 
the 2-dimensional vector {g^{e;x — S2), S2 G 5*2}; the matrices Q^.' (e)s2,si were defined 



13 



after Eq. (l2.16p . Clearly, the whole dependence on x is contained in V^{e;x), whose 
components are "waves" of the same shape (independent of ti,di) originating at the two 
points of 5*2. The shape of these waves is given, far from the source, by the levels of the 
function at the exponent of Eq. fl3.1ip . 

The matrix mtr(e), depending on the details of the interaction v (in particular on di,ti), 
controls the interference of the two waves. Indeed, Eq. (I3.14p describes the squared norm 
of the linear combination mtr(e)^/^l^_(e; x). 

i = 2 (reflected) 

5(2) (e;x) = P(e)o,o - 23i(r_(e,x),Qj'')(e)V^(°)(e;a;))^^ + (Vl(e, x), m,ef (e)r__(e, x))^^ , 

(3.16) 

where 

m,ef(e) = Qj''^(e) ■ (P(e)|5,x5j ■ Q^-''\e) > 0, (3.17) 

with V^^\e;x) denoting the 2-dimensional vector {P{e)x,s2i ^2 € 82}- The first term in 
Eq. fl3.16p is a constant, giving the density per energy and site of "free" electrons of 
energy e in the initial equilibrium state of the reservoir L2. The last term has the same 
significance as for Eq. (I3.14p with mref instead of rritj.. The middle term is new and gives 
account of the overlapping between the "waves" originating in S2 with the free-electron 
states in L2 of energy e. 




Figure 3: Density plot of 5^^\0.3;x) for t2 = 1, di = 1 and: a. ti = 1; b. ti = 1/2; c. ti = 0; 
Legend values: k ■ 10~^, k = 1, 10. 

We start by exploring the "fixed-energy" density distributions, S^^^ and S^'^\ In Figures [3] 
andlH density plots of S^^\0.3;x) and S^'^\0.3;x), respectively, are shown for x in a square of 
side 40 containing 5*2. Thereby, the cases ti = 1 (symmetric), ti = 1/2 (asymmetric) and ti = 
(corresponding to one contact) are represented, in order to exemplify the dependence on ti. 
Figure [S] shows the density plot of S^^\lA; x) in the symmetric case ti = 1 for a. di = 1 and b. 
di = 20; panel c is the plot of 5(^^(1.4; x) in the case of one contact {ti = 0). 
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Figure 4: Density plot of 6^'^\0.3;x) for t2 = I, di = 1 and: a. h = 1; b. h = 1/2; c. ti = 0. 
Legend values: 0.095 + k ■ 0.012, = 1, 10. 




Figure 5: Density plot of S^^\lA;x) for t2 = I and: a. ti = 1, (ii = 1; b. ti = 1, di = 20; c. 
ti = 0. Legend values: 6k ■ 10~^, k = 1, 10. 

The interference patterns are clearly visible in the symmetric case and di = 1 (left panels), 
but the ones for S^^^ and for S^'^^ are drastically different. The number of fringes increases with 
increasing e. In the case of S^^\ the contrast of the fringes decreases as ti decreases, leading, 
when ti = 0, to their complete disappearance (panels c in Figures [S] and H]) . What concerns 
^(2)^ "circular" fringes around S2,2 survive even at ti = (Figure [Sb). 

This can be accounted for using the asymptotic form of V-{e,x), see Eq. (13. lip . Indeed, 
for 5^^\ the interference of the two components g^(e]X — S2), S2 € S2 in the linear combina- 
tion mtr{ey^'^V-{e; x) is controlled by the phase difference \x — S2,i|-ft's,i cos ((/)i — 9s,i) — \x — 
•52,2 1 -f^s, 2 cos (02 — ^5,2), where we denoted with subscripts 1,2 the corresponding functions cal- 
culated for S2,i, S2,2, respectively; taking into account that Ks cos (0 — 6s) is slowly varying with 
(e.g. at e = 0.3, it oscillates between 0.784 and 0.795, and its derivative does not exceed 
0.022), the above difference is constant (modulo 2k7r) on a family of curves similar to hyper- 
bolae with foci S2,i,S2,2- The pattern for S^'^^ comes from the second term in Eq. (I3.16P : as 
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V^^\e;x) is a real vector, this term is a sum of cosines, cos{ipi + \x — S2,j|-ft's,i cos (0j — Og^i)), 
which is constant on curves similar to circles centered at S2,i,S2,2- A similar explanation is 
valid in the situation investigated in [5] and, in fact, the density plots of Figure H] are quite 
similar with those reported in [5]. The number of interference fringes depends on e by the 
monotonicity of the exponent in Eq. fl3.1ip . The fact that the density is higher along diagonals 
of the lattice, which is striking in Figure [5l is due to the fact that the amplitude ipi^y 0) has a 
minimum (maximum) at = 7r/4(resp 0); at e = 1.4, the ratio of \g-\'^ on the diagonal and on 
the axis (at the same distance from the origin) equals ipilA, 0)/'?/'(1.4, vr/4) = 2.264; at e = 0.3, 
this ratio equals only 1.127. This explains, at least for one contact (ti = 0), why the density is 
larger on the diagonals and why the effect is better seen at the larger energy (Figures [3t, Ht). 

One expects a dependence of S^^\e;x) on the distance di between the two points of 5*1, 
as the larger di, the less correlated are the electrons incident at two distant contacts. The 
calculations show that the visibility of the fringes decreases with increasing di (see Figure O 
panel b). However, the interference pattern is conserved even in the limit di — )■ oo. Indeed, 
in this limit P(e)|5'^xS'i becomes diagonal, hence, mtr(e) 
which is by no means diagonal. 



(1,2), 





0.0515 
0.0510 
0.0505 
0.0500 




Figure 6: Local densities on the line x = {i,19),i = 1,...,40: a. Transmitted electrons (/^i = 
1.4); b. Reflected electrons {fi2 = 0.3); c. The (total) local density in the stationary state. 
Dashed, the initial equilibrium density is represented. 

The density distribution in L2 in the stationary state is the sum d{x) = J^'^ 5^^\e; x)de + 
^{2)^g. As position (and number) of the fringes in 5^^\e, x) depends on e, we expect 

the interference pattern of d{x) to have worse visibility. This is confirmed by the computation 
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Figure 7: Graph of j(e) 



of d{x). We fixed, as before, temperature to and /ii = 1.4, /i2 = 0.3. In a density plot of 
d{x), fringes are indeed hardly visible. Indeed, in the intermediate region where fringes are more 

(2) 

pronounced, d{x) is close to the initial equilibrium density of the second reservoir, peq = 0.0492 
and has very small oscillations (of order 10"'^) around this value. In order to better exhibit the 
existence of fringes, we plotted the two terms and their sum on the line x = {i, 19), i = 1, 40 
in Figure [6l 

3.3 The particle current 

The total current flowing from reservoir 1 is given by Eq. f l2.24p . Under the assumptions of the 
previous subsection, in particular Pi = P2 = 00, fii = 1.4, fi2 = 0.3, it is easy to see that 

1.4 

Jl, = j j(e)de; = ^-kIi^ [m,,{e)P{e)\s,^s,] , (3.18) 
0.3 

with mtr(e) the same matrix as in Eq. (13.151) . We report the calculation for the symmetric case 
^1 = ^2 = 1 and di = 1, d2 = 20: Jl^ = 0.2416. It is interesting to exhibit the energy resolution 
of the current, as well. The plot of j(e), represented in Figure [TJ shows oscillations, indicating 
that tuning the energy may result in higher conductivities. The dotted line in Figure [7] is the 
plot of 2jo(e), where jo{e) corresponds to one contact (ti = 0, ^2 = 1, or viceversa). The latter 
plot, which is what Ohm's law would predict, has no structure and definitely exceeds j'(e) over 
the whole domain. 

Interference effects are clearly seen in the distribution of the local currents in a neighborhood 
of the contacts. We computed, using Eq. fl2.20p . the local currents in the stationary state in 
the second reservoir in the same 40 x 40 square around 5*2, for the same setting and parameter 
values as before. For {x,?/} a bond of nearest neighbors in L2, one has, after integration over 
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Figure 8: The local "transmitted" current: a. at energy e = 1.4; b. at energy e = 0.3; c. 
integrated on the energy range (0, 1.4), across the vertical bonds {{i, 19), (i, 20)}. Dashed, the 
local density in points (i, 19) is plotted. 
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Figure 9: The local currents in the stationary state across the line of bonds {{i, 19), {i, 20)}, i 
1, ...,40. Dotted is the local density (shifted downwards with 0.048) at sites {i, 19). 



the energy shells. 



0.3 



1.4 



lx,y 



Jx,y 



;e)de + j 



Jx,y 



e)de, 



where: 



jS(e) 



-%V^{e, x), mtr(e)Vl(e, y))^2, 



- ^^(K_(e, x), mref (e)F_(e, y))^ 



(3.19) 



Figure [8] is a plot of the local currents of transmitted electrons along nearest-neighbor 
bonds crossing a horizontal line in the region of fringes: a. of fixed energy e = 1.4, i.e. 
(1.4); b. of fixed energy e = 0.3, (0.3); c. integrated over the whole energy range. 



.(tr) 
I.e. ]x,y 



1.4 



/ (e)de. We took {x, ?/} = {(i, 19), (i, 20)}, % = 1,...,40. We represented also 



the local density on that line, in order to demonstrate the correlation of the two quantities: 
currents are larger on bonds starting from high-density sites. The "reflected" local currents 
jx!y^ across the same line show the same kind of correlation with the local density, the major 
difference being that the currents are negative, i.e. electrons enter the neighborhood of 5*2. 

Finally, we plotted in Figure [9] the total (transmitted + reflected) local currents in the 
stationary state, jx,y, on the same line of bonds {x,y} = {{i, 19), (i,20)}, i = I, ...,40, along 
with a copy of the plot of the local density at x (Figure [6t). As the variation of jx,y is an order 
of magnitude larger than that of jx,f\ the picture is quite similar to that of jx^J (Figure [8t), 
and the same correlation with the total local density is observed. 



4 Conclusion 



We have performed a detailed study of the stationary state, for the model under consideration, 
by calculating the expectations of various local observables, namely, the number of particles in 
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a lattice site and the particle current along a nearest-neighbor bond (in a particular geometric 
arrangement and tunneling constants). These expectations exhibit a peculiar dependence on 
the position of the lattice site/bond relative to the two contacts, putting into evidence inter- 
ference patterns consisting of fringes of high local density, respectively local current. We chose 
to calculate this space dependence in the less populated of the two systems. Both the density 
and current profiles are sums of contributions from transmitted and reflected particles (corre- 
sponding to one of the systems being unpopulated). The two contributions yield qualitatively 
different patterns: the fringes of the transmitted particles are similar to hyperbolae with foci 
in the contacts and they disappear as one of the contacts is suppressed, while the fringes of the 
reflected particles are similar to circles centered in the contacts and fringes are present even in 
the case of one contact. We explained this peculiarity, as well as the dependence of the pattern 
on the various parameters. The space dependence of the sum of the two contributions turns 
out, in our case, to be dominated by the transmitted particles: hyperbolic fringes of higher 
density along which particles flow from the contacts to infinity. 

New interference effects are to be expected when having, instead of direct tunneling, some 
finite intermediate quantum system, possibly with a Coulomb repulsion included, in which case 
the scattering has a resonant structure. For small tunneling constants, resonances come close to 
the energy levels of the intermediate quantum system, with the effect that only certain energy 
channels are open. In this case, the fixed-energy calculations performed above become relevant. 
We propose to follow the subject in another publication. 
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